Esta práctica cubre de forma transversal la asignatura.
Las Prácticas 1 y 2 de la asignatura se plantean de una forma conjunta de modo que la Práctica 2 será continuación de la 1.
El objetivo global de las dos prácticas consiste en seleccionar uno o varios juegos de datos, realizar las tareas de preparación y análisis exploratorio con el objetivo de disponer de datos listos para aplicar algoritmos de clustering, asociación y clasificación.
Las competencias que se trabajan en esta prueba son:
La correcta asimilación de todos los aspectos trabajados durante el semestre.
En esta práctica abordamos un caso real de miner??a de datos donde tenemos que poner en juego todos los conceptos trabajados. Hay que trabajar todo el ciclo de vida del proyecto. Desde el objetivo del proyecto hasta la implementación del conocimiento encontrado pasando por la preparación, limpieza de los datos, conocimiento de los datos, generación del modelo, interpretación y evaluación.
Material docente proporcionado por la UOC.
Ejercicios prácticos
Para todas las PEC es necesario documentar en cada apartado del ejercicio práctico que se ha hecho y como se ha hecho.
El formato de entrega es: usernameestudiante-PECn.html/doc/docx/odt/pdf
Fecha de entrega: 15/01/2020
Se debe entregar la PEC en el buzón de entregas del aula
A menudo es inevitable, al producir una obra multimedia, hacer uso de recursos creados por terceras personas. Es por lo tanto comprensible hacerlo en el marco de una práctica de los estudios de Informática, Multimedia y Telecomunicación de la UOC, siempre y cuando esto se documente claramente y no suponga plagio en la práctica.
Por lo tanto, al presentar una práctica que haga uso de recursos ajenos, se debe presentar junto con ella un documento en que se detallen todos ellos, especificando el nombre de cada recurso, su autor, el lugar donde se obtuvo y su estatus legal: si la obra esta protegida por el copyright o se acoge a alguna otra licencia de uso (Creative Commons, licencia GNU, GPL …). El estudiante deberá asegurarse de que la licencia no impide espec??ficamente su uso en el marco de la práctica. En caso de no encontrar la información correspondiente tendrá que asumir que la obra esta protegida por copyright.
Deberéis, además, adjuntar los ficheros originales cuando las obras utilizadas sean digitales, y su código fuente si corresponde.
Como continuación del estudio iniciado en la práctica 1, procedemos en esta práctica 2 a aplicar modelos anal??ticos sobe el juego de datos seleccionado y preparado en la práctica anterior.
De este modo se pide al estudiante que complete los siguientes pasos:
Aplicar un modelo de generación de reglas a partir de Reglas de asociacion (Hay que usar el algoritmo arules).
Aplicar un modelo no supervisado y basado en el concepto de distancia, sobre el juego de datos (Hay que usar el algoritmo kmeans).
Aplica de nuevo el modelo anterior, pero usando una métrica distinta y compara los resultados (Hay que usar un algoritmo de clustering diferente al kmeans)
Aplicar un modelo supervisado sobre el juego de datos sin haber aplicado previamente PCA/SVD (Hay que usar un árbol de decisión).
Aplicar un modelo supervisado sobre el juego de datos habiendo aplicado previamente PCA/SVD (Hay que usar un árbol de decisión tras la aplicación de un PCA/SVD).
¿Ha habido mejora en capacidad predictiva, tras aplicar PCA/SVD? ¿A qué crees que es debido? (Hay que comparar los resultados de los puntos 4 y 5).
Vamos a comenzar realizando la carga de las librerías y funciones necesarias para ejecutar el ejercicio.
libraries <- c("caret", "arules", "cluster", "fpc","tibble", "factoextra", "C50", "ggbiplot", "dplyr", "cluster")
check.libraries <- is.element(libraries, installed.packages()[, 1]) == FALSE
libraries.to.install <- libraries[check.libraries]
if (length(libraries.to.install != 0)) {
install.packages(libraries.to.install)
}
success <- sapply(libraries, require, quietly = FALSE, character.only = TRUE)
if(length(success) != length(libraries)) {
stop("A package failed to return a success in require() function.")
}normalize <- function(x) {
return ((x - min(x)) / (max(x) - min(x)))
}De la resolución de la práctica anterior, se han exportado los datos finales y vueltos a importar dentro del contexto de este proyecto.
data.raw <- read.csv('dengue_data.csv',stringsAsFactors = FALSE)
data.raw$cases <- as.factor(data.raw$cases)
summary(data.raw)## X ID city year
## Min. : 1.0 Min. : 1.0 Length:1456 Min. :1990
## 1st Qu.: 364.8 1st Qu.: 364.8 Class :character 1st Qu.:1997
## Median : 728.5 Median : 728.5 Mode :character Median :2002
## Mean : 728.5 Mean : 728.5 Mean :2001
## 3rd Qu.:1092.2 3rd Qu.:1092.2 3rd Qu.:2005
## Max. :1456.0 Max. :1456.0 Max. :2010
## weekofyear trimester ndvi precipitation
## Min. : 1.00 Length:1456 Min. :0.0000 Min. :0.0000
## 1st Qu.:13.75 Class :character 1st Qu.:0.3239 1st Qu.:0.0255
## Median :26.50 Mode :character Median :0.4085 Median :0.0991
## Mean :26.50 Mean :0.4356 Mean :0.1172
## 3rd Qu.:39.25 3rd Qu.:0.5353 3rd Qu.:0.1793
## Max. :53.00 Max. :1.0000 Max. :1.0000
## air_temp_mean diur_temp_range humidity total_cases
## Min. :0.0000 Min. :0.00000 Min. :0.0000 Min. : 0.00
## 1st Qu.:0.4005 1st Qu.:0.06621 1st Qu.:0.4397 1st Qu.: 5.00
## Median :0.5323 Median :0.10273 Median :0.6116 Median : 12.00
## Mean :0.5375 Mean :0.24174 Mean :0.5752 Mean : 24.68
## 3rd Qu.:0.6864 3rd Qu.:0.42235 3rd Qu.:0.7157 3rd Qu.: 28.00
## Max. :1.0000 Max. :1.00000 Max. :1.0000 Max. :461.00
## cases
## low :1431
## moderate: 15
## severe : 10
##
##
## str(data.raw)## 'data.frame': 1456 obs. of 13 variables:
## $ X : int 1 2 3 4 5 6 7 8 9 10 ...
## $ ID : int 1 2 3 4 5 6 7 8 9 10 ...
## $ city : chr "sj" "sj" "sj" "sj" ...
## $ year : int 1990 1990 1990 1990 1990 1990 1990 1990 1990 1990 ...
## $ weekofyear : int 18 19 20 21 22 23 24 25 26 27 ...
## $ trimester : chr "T2" "T2" "T2" "T2" ...
## $ ndvi : num 0.408 0.419 0.379 0.506 0.556 ...
## $ precipitation : num 0.0318 0.0584 0.0884 0.0393 0.0193 ...
## $ air_temp_mean : num 0.388 0.473 0.548 0.575 0.646 ...
## $ diur_temp_range: num 0.0867 0.0691 0.0643 0.073 0.113 ...
## $ humidity : num 0.263 0.418 0.587 0.567 0.628 ...
## $ total_cases : int 4 5 4 3 6 2 4 5 10 6 ...
## $ cases : Factor w/ 3 levels "low","moderate",..: 1 1 1 1 1 1 1 1 1 1 ...Dado que estos datos provienen de un ejercicio en el que se ha elaborado una amplia y exhaustiva tarea de limpieza y normalización de datos, vamos a saltarnos este paso y pasar directamente a la elaboración de los modelos.
Sin embargo, existe un problema que debe ser solventado y que ha sido descubierto durante la fase preliminar, y no documentada, de la elaboración del ejercicio. El problema surgió tras explorar la generación de reglas mediante arules cuando el algoritmo no generaba regla alguna, o lo que es lo mismo, todas las reglas estaban orientadas a la misma clase resultado.
Tras múltiples intentos y evaluaciones se observó que la variable resultante cases tiene una distribución poco equilibrada y no apta para la generación de reglas. Como era de esperar, este problema no estaba presente en los modelos no supervisados, aunque todos los modelos supervisamos daban como resultado una precisión superior al 98%, algo poco común.
Aunque tiene sentido agrupar los casos de dengue de forma que a priori nos parezca razonable, la realidad es que desde un punto de vista estadístico no lo es. Si observamos, más del 95% de los casos de dengue caen bajo la categoría de low o bajos. Esto hace que cualquier modelo generado a partir de estos datos sea poco o nada eficaz. El modelo será muy eficiente prediciendo la baja indicencia de dengue pero nefasto en la predicción de casos moderados o de alta incidencia de la enfermedad.
Para resolver este problema de distribución de las clases resultantes, se ha optado por la redistribución de esta variable siguiendo otros criterios más objetivos. Para ello se ha hecho uso de la media como línea de corte; es decir, que todas aquellas semanas donde se hayan detectado hasta 25 casos de dengue se considerará una semana normal, mientras que si superamos los 25 casos lo catalogaremos como pandemic.
summary(data.raw$total_cases)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 5.00 12.00 24.68 28.00 461.00Tras realizar la nueva discretización de los datos mediante el método fixed, objetemos una distribución de 1036 casos normales y 420 pandémicos. Esta nueva distribución puede, sin duda alguna, generar modelos más fiables en su predicción.
outcome <- data.frame('outcome' = discretize(data.raw$total_cases,method="fixed", breaks = c(0,25, 500), labels=c('normal', 'pandemic')))
data.raw <- add_column(data.raw, outcome, .after = 13)
data.raw$cases <- NULL
summary(data.raw)## X ID city year
## Min. : 1.0 Min. : 1.0 Length:1456 Min. :1990
## 1st Qu.: 364.8 1st Qu.: 364.8 Class :character 1st Qu.:1997
## Median : 728.5 Median : 728.5 Mode :character Median :2002
## Mean : 728.5 Mean : 728.5 Mean :2001
## 3rd Qu.:1092.2 3rd Qu.:1092.2 3rd Qu.:2005
## Max. :1456.0 Max. :1456.0 Max. :2010
## weekofyear trimester ndvi precipitation
## Min. : 1.00 Length:1456 Min. :0.0000 Min. :0.0000
## 1st Qu.:13.75 Class :character 1st Qu.:0.3239 1st Qu.:0.0255
## Median :26.50 Mode :character Median :0.4085 Median :0.0991
## Mean :26.50 Mean :0.4356 Mean :0.1172
## 3rd Qu.:39.25 3rd Qu.:0.5353 3rd Qu.:0.1793
## Max. :53.00 Max. :1.0000 Max. :1.0000
## air_temp_mean diur_temp_range humidity total_cases
## Min. :0.0000 Min. :0.00000 Min. :0.0000 Min. : 0.00
## 1st Qu.:0.4005 1st Qu.:0.06621 1st Qu.:0.4397 1st Qu.: 5.00
## Median :0.5323 Median :0.10273 Median :0.6116 Median : 12.00
## Mean :0.5375 Mean :0.24174 Mean :0.5752 Mean : 24.68
## 3rd Qu.:0.6864 3rd Qu.:0.42235 3rd Qu.:0.7157 3rd Qu.: 28.00
## Max. :1.0000 Max. :1.00000 Max. :1.0000 Max. :461.00
## outcome
## normal :1036
## pandemic: 420
##
##
##
## Normalizamos la columna weekofyears para poder hacer uso de esta variable en un futuro como variable discreata.
weekofyear_norm <- as.data.frame(apply(data.raw[which( colnames(data.raw)=="weekofyear" )], 2, normalize))
data.raw$weekofyear <- NULL
data.raw <- add_column(data.raw, weekofyear_norm, .after = 4)
summary(data.raw)## X ID city year
## Min. : 1.0 Min. : 1.0 Length:1456 Min. :1990
## 1st Qu.: 364.8 1st Qu.: 364.8 Class :character 1st Qu.:1997
## Median : 728.5 Median : 728.5 Mode :character Median :2002
## Mean : 728.5 Mean : 728.5 Mean :2001
## 3rd Qu.:1092.2 3rd Qu.:1092.2 3rd Qu.:2005
## Max. :1456.0 Max. :1456.0 Max. :2010
## weekofyear trimester ndvi precipitation
## Min. :0.0000 Length:1456 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.2452 Class :character 1st Qu.:0.3239 1st Qu.:0.0255
## Median :0.4904 Mode :character Median :0.4085 Median :0.0991
## Mean :0.4905 Mean :0.4356 Mean :0.1172
## 3rd Qu.:0.7356 3rd Qu.:0.5353 3rd Qu.:0.1793
## Max. :1.0000 Max. :1.0000 Max. :1.0000
## air_temp_mean diur_temp_range humidity total_cases
## Min. :0.0000 Min. :0.00000 Min. :0.0000 Min. : 0.00
## 1st Qu.:0.4005 1st Qu.:0.06621 1st Qu.:0.4397 1st Qu.: 5.00
## Median :0.5323 Median :0.10273 Median :0.6116 Median : 12.00
## Mean :0.5375 Mean :0.24174 Mean :0.5752 Mean : 24.68
## 3rd Qu.:0.6864 3rd Qu.:0.42235 3rd Qu.:0.7157 3rd Qu.: 28.00
## Max. :1.0000 Max. :1.00000 Max. :1.0000 Max. :461.00
## outcome
## normal :1036
## pandemic: 420
##
##
##
## Como paso extra, vamos a buscar correlaciones en los datos para evaluar si es posible deshacernos de alguna variable. Los datos indican que no existe correlación entre los datos ya que ninguno es próximo a -1 o 1.
cor(data.raw[,c(5,7,8,9,10,11)])## weekofyear ndvi precipitation air_temp_mean
## weekofyear 1.00000000 0.04094926 0.11662773 0.42134406
## ndvi 0.04094926 1.00000000 0.17711009 -0.30325048
## precipitation 0.11662773 0.17711009 1.00000000 -0.01721905
## air_temp_mean 0.42134406 -0.30325048 -0.01721905 1.00000000
## diur_temp_range 0.07712561 0.67396605 0.20132179 -0.28037539
## humidity 0.34718816 0.07904915 0.45060559 0.50838016
## diur_temp_range humidity
## weekofyear 0.07712561 0.34718816
## ndvi 0.67396605 0.07904915
## precipitation 0.20132179 0.45060559
## air_temp_mean -0.28037539 0.50838016
## diur_temp_range 1.00000000 0.01294761
## humidity 0.01294761 1.00000000Para generar las reglas vamos a trabajar con dos subconjuntos de datos. En el primero utilizaremos la variable trimestre, que como ya sabemos, es una indicadora de las estaciones del año, y en el segundo trabajaremos sólo con datos climatológicos. La razón por la que hemos decidido utilizar dos subconjuntos de datos es porque creemos que el dato temporal es, en cierto modo, una función de los datos climatológicos y morfológicos recogidos.
Queremos ver la influencia directa que la variable trimestre pueda tener sobre las reglas en contraste con los factores climatológicos y morfológicos por si solos, y lo más importante, estudiar la eficacia de los modelos generados de ambos conjuntos.
data.arules.seasons <- data.raw[, -c(1, 2, 3, 4, 5, 12)]
summary(data.arules.seasons)## trimester ndvi precipitation air_temp_mean
## Length:1456 Min. :0.0000 Min. :0.0000 Min. :0.0000
## Class :character 1st Qu.:0.3239 1st Qu.:0.0255 1st Qu.:0.4005
## Mode :character Median :0.4085 Median :0.0991 Median :0.5323
## Mean :0.4356 Mean :0.1172 Mean :0.5375
## 3rd Qu.:0.5353 3rd Qu.:0.1793 3rd Qu.:0.6864
## Max. :1.0000 Max. :1.0000 Max. :1.0000
## diur_temp_range humidity outcome
## Min. :0.00000 Min. :0.0000 normal :1036
## 1st Qu.:0.06621 1st Qu.:0.4397 pandemic: 420
## Median :0.10273 Median :0.6116
## Mean :0.24174 Mean :0.5752
## 3rd Qu.:0.42235 3rd Qu.:0.7157
## Max. :1.00000 Max. :1.0000A continuación, discretizamos las variables numéricas.
data.arules.seasons$ndvi <- discretize(data.arules.seasons$ndvi, method="interval", labels=c('water', 'cement', 'vegetal'))
data.arules.seasons$precipitation <- discretize(data.arules.seasons$precipitation, method="interval", labels=c('dry', 'normal', 'storm'))
data.arules.seasons$air_temp_mean <- discretize(data.arules.seasons$air_temp_mean, method="interval", labels=c('low', 'medium', 'high'))
data.arules.seasons$diur_temp_range <- discretize(data.arules.seasons$diur_temp_range, method="interval", labels=c('low', 'medium', 'high'))
data.arules.seasons$humidity <- discretize(data.arules.seasons$humidity, method="interval", labels=c('dry', 'normal', 'wet'))summary(data.arules.seasons)## trimester ndvi precipitation air_temp_mean diur_temp_range
## Length:1456 water :399 dry :1398 low :184 low :1000
## Class :character cement :919 normal: 55 medium:872 medium: 332
## Mode :character vegetal:138 storm : 3 high :400 high : 124
## humidity outcome
## dry :158 normal :1036
## normal:738 pandemic: 420
## wet :560Creamos los conjuntos de entramiento y test.
set.seed(777)
train <- createDataPartition(data.arules.seasons$outcome, p = 0.66, list = FALSE)
data.arules.seasons.train <- data.arules.seasons[train,]
data.arules.seasons.test <- data.arules.seasons[-train,]
data.arules.seasons.train.x <- data.arules.seasons.train[,1:6]
data.arules.seasons.train.y <- data.arules.seasons.train[,7]
data.arules.seasons.test.x <- data.arules.seasons.test[,1:6]
data.arules.seasons.test.y <- data.arules.seasons.test[,7]Y estudiamos los resultados del modelo generado.
data.arules.seasons.model <- C5.0(x = data.arules.seasons.train.x, y = data.arules.seasons.train.y, rules=TRUE)
summary(data.arules.seasons.model)##
## Call:
## C5.0.default(x = data.arules.seasons.train.x, y =
## data.arules.seasons.train.y, rules = TRUE)
##
##
## C5.0 [Release 2.07 GPL Edition] Mon Jun 8 20:53:17 2020
## -------------------------------
##
## Class specified by attribute `outcome'
##
## Read 962 cases (7 attributes) from undefined.data
##
## Rules:
##
## Rule 1: (290/19, lift 1.3)
## diur_temp_range in {medium, high}
## -> class normal [0.932]
##
## Rule 2: (482/78, lift 1.2)
## trimester in {T2, T1}
## -> class normal [0.837]
##
## Rule 3: (322/132, lift 2.0)
## trimester in {T3, T4}
## diur_temp_range = low
## -> class pandemic [0.590]
##
## Default class: normal
##
##
## Evaluation on training data (962 cases):
##
## Rules
## ----------------
## No Errors
##
## 3 220(22.9%) <<
##
##
## (a) (b) <-classified as
## ---- ----
## 552 132 (a): class normal
## 88 190 (b): class pandemic
##
##
## Attribute usage:
##
## 83.58% trimester
## 63.62% diur_temp_range
##
##
## Time: 0.0 secsLas reglas encontradas indican que:
Rango en temperatura diurna (medio, alto) -> casos normales
Estaciones de invierno y primavera -> casos normales
Estaciones verano y otoño + Rango en temperatura diurna (bajo) -> casos pandémicos.
Según estas reglas, existen dos factores que influencian la proliferación de los casos de dengue. Estos son las estaciones del año y la temperatura media diurna. Los casos aumentan durante el verano y el otoño (estación de huracanes para el área de San juan) y con un rango de temperatura media bajo; esto es, con pocas o ninguna bajada de temperatura durante la tarde. Los casos normales de dengue ocurren bajo condiciones totalmente opuestas. Podemos decir que al mosquito le gusta el calor y el agua.
A pesar de ello, no vemos que el modelo pueda describir muy bien las condiciones exactas bajo las que prolifera el mosquito ya que sólo describe el 50% de los casos pandémicos. También podemos observar que el modelo deja la predicción en manos de la estación del año ignorando el resto de las variables.
data.arules.seasons.model.chart <- C5.0(x = data.arules.seasons.train.x, y = data.arules.seasons.train.y)
plot(data.arules.seasons.model.chart, subtree = 1)
Vamos a proceder ahora con el segundo conjunto de datos donde ignoramos la estación del año.
data.arules.no_seasons <- data.arules.seasons[,-1]
summary(data.arules.no_seasons)## ndvi precipitation air_temp_mean diur_temp_range humidity
## water :399 dry :1398 low :184 low :1000 dry :158
## cement :919 normal: 55 medium:872 medium: 332 normal:738
## vegetal:138 storm : 3 high :400 high : 124 wet :560
## outcome
## normal :1036
## pandemic: 420
## Creamos los conjuntos de entramiento y test.
set.seed(777)
train <- createDataPartition(data.arules.no_seasons$outcome, p = 0.66, list = FALSE)
data.arules.no_seasons.train <- data.arules.no_seasons[train,]
data.arules.no_seasons.test <- data.arules.no_seasons[-train,]
data.arules.no_seasons.train.x <- data.arules.no_seasons.train[,1:5]
data.arules.no_seasons.train.y <- data.arules.no_seasons.train[,6]
data.arules.no_seasons.test.x <- data.arules.no_seasons.test[,1:5]
data.arules.no_seasons.test.y <- data.arules.no_seasons.test[,6]Y estudiamos los resultados del modelo generado.
data.arules.no_seasons.model <- C5.0(x = data.arules.no_seasons.train.x, y = data.arules.no_seasons.train.y, rules=TRUE)
summary(data.arules.no_seasons.model)##
## Call:
## C5.0.default(x = data.arules.no_seasons.train.x, y
## = data.arules.no_seasons.train.y, rules = TRUE)
##
##
## C5.0 [Release 2.07 GPL Edition] Mon Jun 8 20:53:17 2020
## -------------------------------
##
## Class specified by attribute `outcome'
##
## Read 962 cases (6 attributes) from undefined.data
##
## Rules:
##
## Rule 1: (290/19, lift 1.3)
## diur_temp_range in {medium, high}
## -> class normal [0.932]
##
## Rule 2: (693/147, lift 1.1)
## air_temp_mean in {low, medium}
## -> class normal [0.787]
##
## Rule 3: (254/124, lift 1.8)
## air_temp_mean = high
## diur_temp_range = low
## -> class pandemic [0.512]
##
## Default class: normal
##
##
## Evaluation on training data (962 cases):
##
## Rules
## ----------------
## No Errors
##
## 3 272(28.3%) <<
##
##
## (a) (b) <-classified as
## ---- ----
## 560 124 (a): class normal
## 148 130 (b): class pandemic
##
##
## Attribute usage:
##
## 98.44% air_temp_mean
## 56.55% diur_temp_range
##
##
## Time: 0.0 secsAquí obtenemos las siguientes reglas:
Rango en temperatura diurna (medio, alto) -> casos normales
Temperatura aire promedio (bajo, medio) -> casos normales
Rango en temperatura diurna (bajo) + Temperatura aire promedio (alto) -> casos pandémicos
Para este conjunto de datos, el árbol resultante es más diverso al proveer reglas que incluyen otras variables, podemos decir que “se explica mejor” que el modelo anterior. Aquí podemos ver que los casos de dengue aumentan con las combinaciones de altas temperaturas promedios diarias y constantes durante la tarde y sobre todo en zonas urbanas (ndvi ~ 0).
El factor temperatura ha sido constante en ambos estudios, extendido los resultados a la temperatura media en este segundo conjunto de datos donde además se ha sumado el factor ndvi en forma de zonas con cemento o urbanizadas, que es precisamente donde el agua es más propensa a estancarse. A esta informacion, y según el modelo anterior, podemos añadirle el hecho de que las estaciones de verano y otoño son las mas propensas para desarrollo del mosquito.
data.arules.no_seasons.model.chart <- C5.0(x = data.arules.no_seasons.train.x, y = data.arules.no_seasons.train.y)
plot(data.arules.no_seasons.model.chart, subtree = 1)
Vamos a continuación a comparar ambos modelos en términos de eficiencia predictiva de los datos de prueba.
data.arules.seasons.predicted_model <- predict(data.arules.seasons.model, data.arules.seasons.test.x, type="class")
print(sprintf("La precisión del árbol con estaciones es: %.4f %%",100*sum(data.arules.seasons.predicted_model == data.arules.seasons.test.y) / length(data.arules.seasons.predicted_model)))## [1] "La precisión del árbol con estaciones es: 78.9474 %"data.arules.no_seasons.predicted_model <- predict(data.arules.no_seasons.model, data.arules.no_seasons.test.x, type="class")
print(sprintf("La precisión del árbol sin estaciones es: %.4f %%",100*sum(data.arules.no_seasons.predicted_model == data.arules.no_seasons.test.y) / length(data.arules.no_seasons.predicted_model)))## [1] "La precisión del árbol sin estaciones es: 71.2551 %"El árbol deducido del conjunto de datos que incluye las estaciones del año ha devengado un mayor porcentaje en la precisión de los resultados en contraste con aquel que utiliza datos puramente climáticos y morfológicos. Por el contrario, la definición de reglas es más rica y explicita al omitir la variable estacional.
Para realizar el modelo no supervisado vamos a repetir el procedimiento anterior y hacer uso de un conjunto con variable temporal y otro sin ella. Utilizaremos la variable normalizada weekofyear_morm. El propósito de este ejercicio es encontrar el número de clústeres óptimo para nuestro modelo mediante el método kmeans.
data.kmeans.seasons <- data.raw[, c(5,7,8,9,10,11,13)]
summary(data.kmeans.seasons)## weekofyear ndvi precipitation air_temp_mean
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.2452 1st Qu.:0.3239 1st Qu.:0.0255 1st Qu.:0.4005
## Median :0.4904 Median :0.4085 Median :0.0991 Median :0.5323
## Mean :0.4905 Mean :0.4356 Mean :0.1172 Mean :0.5375
## 3rd Qu.:0.7356 3rd Qu.:0.5353 3rd Qu.:0.1793 3rd Qu.:0.6864
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.0000
## diur_temp_range humidity outcome
## Min. :0.00000 Min. :0.0000 normal :1036
## 1st Qu.:0.06621 1st Qu.:0.4397 pandemic: 420
## Median :0.10273 Median :0.6116
## Mean :0.24174 Mean :0.5752
## 3rd Qu.:0.42235 3rd Qu.:0.7157
## Max. :1.00000 Max. :1.0000Para el cálculo de clusters óptimo, realizamos una iteración en la que asumimos un total de 2 hasta de 10 clusters. Con cada uno de estos valores utilizamos el algoritmo kmeans para clasificar las muestras y el método silhouette para estudiar la precisión en la agrupación de los datos.
x <- data.kmeans.seasons[,1:6]
d <- daisy(x)
resultados <- rep(0, 7)
cluster <- data.frame(c(2,3,4,5,6,7))
for (i in c(2,3,4,5,6,7))
{
fit <- kmeans(x, i)
y_cluster <- fit$cluster
sk <- silhouette(y_cluster, d)
resultados[i] <- mean(sk[,3])
}
plot(2:7,resultados[2:7],type="o",col="blue",pch=0,xlab="Número de clusters",ylab="Silueta")
La gráfica revela una agrupación más definida para un modelo con 4 clústeres. El método elbow, debajo, recomienda el uso de 3 clústeres.
resultados <- rep(0, 7)
for (i in c(2,3,4,5,6,7))
{
fit <- kmeans(x, i)
resultados[i] <- fit$tot.withinss
}
plot(2:7,resultados[2:7],type="o",col="blue",pch=0,xlab="Número de clusters",ylab="tot.tot.withinss")
Si procedemos con la función kmeansruns bajos los criterios de la silueta media y Calinski-Harabasz obtenemos que para ambos casos el número óptimo de clusters es de 4.
fit_ch <- kmeansruns(x, krange = 1:10, criterion = "ch")
fit_asw <- kmeansruns(x, krange = 1:10, criterion = "asw")fit_ch$bestk## [1] 4fit_asw$bestk## [1] 4plot(1:10,fit_ch$crit,type="o",col="blue",pch=0,xlab="Número de clústers",ylab="Criterio Calinski-Harabasz")
plot(1:10,fit_asw$crit,type="o",col="blue",pch=0,xlab="Número de clústers",ylab="Criterio silueta media")
Vamos ahora a realizar las gráficas de los modelos de 2, 3 y 4 clústeres para estudiar la agrupación de los datos en relación con el número de clústeres de forma visual.
Para el clusplot de dos clústeres, que es como tenemos agrupados los datos, podemos ver que aunque si existen dos grupos de datos claramente definidos, las elipses van en dirección contraria a la agrupación natural de los datos. Los mismo ocurre para los gráficos de 3 y 4 clústeres donde las elipses no definen claramente los núcleos de datos.
data.kmeans.2_clusters <- kmeans(x, 2)
clusplot(x, data.kmeans.2_clusters$cluster, color=TRUE, shade=TRUE, labels=2, lines=0)
data.kmeans.3_clusters <- kmeans(x, 3)
clusplot(x, data.kmeans.3_clusters$cluster, color=TRUE, shade=TRUE, labels=3, lines=0)
data.kmeans.4_clusters <- kmeans(x, 4)
clusplot(x, data.kmeans.4_clusters$cluster, color=TRUE, shade=TRUE, labels=3, lines=0)
Al contrastar los datos obtenidos en el análisis de reglas de semana del año y temperatura media podemos ver claramente la polarización de los datos, la cual se repite en los gráficos para 3 y 4 clústeres.
par(mfrow=c(2,2))
plot(data.kmeans.seasons[c(1,4)], col=data.kmeans.2_clusters$cluster)
plot(data.kmeans.seasons[c(1,4)], col=data.kmeans.3_clusters$cluster)
plot(data.kmeans.seasons[c(1,4)], col=data.kmeans.4_clusters$cluster)
El patrón se repite para la variable que define la vegetación de la zona, otra de las reglas relevantes obtenidas.
par(mfrow=c(2,2))
plot(data.kmeans.seasons[c(1,2)], col=data.kmeans.2_clusters$cluster)
plot(data.kmeans.seasons[c(1,2)], col=data.kmeans.3_clusters$cluster)
plot(data.kmeans.seasons[c(1,2)], col=data.kmeans.4_clusters$cluster)
data.kmeans.2_clusters.table <- table(data.kmeans.2_clusters$cluster,data.kmeans.seasons$outcome)
data.kmeans.2_clusters.table##
## normal pandemic
## 1 615 98
## 2 421 322data.kmeans.2_clusters.table.precision = 100*(max(data.kmeans.2_clusters.table[1], data.kmeans.2_clusters.table[3])+ max(data.kmeans.2_clusters.table[2], data.kmeans.2_clusters.table[4]))/(sum(data.kmeans.2_clusters.table))
data.kmeans.2_clusters.table.precision## [1] 71.15385data.kmeans.3_clusters.table <- table(data.kmeans.3_clusters$cluster,data.kmeans.seasons$outcome)
data.kmeans.3_clusters.table##
## normal pandemic
## 1 272 301
## 2 476 98
## 3 288 21data.kmeans.3_clusters.table.precision = 100*(max(data.kmeans.3_clusters.table[1], data.kmeans.3_clusters.table[4]) + max(data.kmeans.3_clusters.table[2], data.kmeans.3_clusters.table[5]) + max(data.kmeans.3_clusters.table[3], data.kmeans.3_clusters.table[6]))/(sum(data.kmeans.3_clusters.table))
data.kmeans.3_clusters.table.precision## [1] 73.1456data.kmeans.4_clusters.table <- table(data.kmeans.4_clusters$cluster,data.kmeans.seasons$outcome)
data.kmeans.4_clusters.table##
## normal pandemic
## 1 278 21
## 2 62 124
## 3 227 178
## 4 469 97data.kmeans.4_clusters.table.precision = 100*(max(data.kmeans.4_clusters.table[1], data.kmeans.4_clusters.table[5])+max(data.kmeans.4_clusters.table[2], data.kmeans.4_clusters.table[6])+max(data.kmeans.4_clusters.table[3], data.kmeans.4_clusters.table[7])+max(data.kmeans.4_clusters.table[4], data.kmeans.4_clusters.table[8]))/(sum(data.kmeans.4_clusters.table))
data.kmeans.4_clusters.table.precision## [1] 75.41209Al calcular la precisión de los tres modelos, aquel que hace uso de los 4 clústeres es el que mayor precisión presenta con un 75.41%. Este valor supera en un 4% y un 2% a los modelos de 2 y 3 clústeres respectivamente.
Repitamos ahora el procedimiento omitiendo la variable temporal de semana del año.
data.kmeans.no_seasons <- data.raw[, c(7,8,9,10,11,13)]
summary(data.kmeans.no_seasons)## ndvi precipitation air_temp_mean diur_temp_range
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.00000
## 1st Qu.:0.3239 1st Qu.:0.0255 1st Qu.:0.4005 1st Qu.:0.06621
## Median :0.4085 Median :0.0991 Median :0.5323 Median :0.10273
## Mean :0.4356 Mean :0.1172 Mean :0.5375 Mean :0.24174
## 3rd Qu.:0.5353 3rd Qu.:0.1793 3rd Qu.:0.6864 3rd Qu.:0.42235
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.00000
## humidity outcome
## Min. :0.0000 normal :1036
## 1st Qu.:0.4397 pandemic: 420
## Median :0.6116
## Mean :0.5752
## 3rd Qu.:0.7157
## Max. :1.0000Inicialmente observamos que la recomendación es de 2 clústeres.
x <- data.kmeans.no_seasons[,1:5]
d <- daisy(x)
resultados <- rep(0, 7)
cluster <- data.frame(c(2,3,4,5,6,7))
for (i in c(2,3,4,5,6,7))
{
fit <- kmeans(x, i)
y_cluster <- fit$cluster
sk <- silhouette(y_cluster, d)
resultados[i] <- mean(sk[,3])
}
plot(2:7,resultados[2:7],type="o",col="blue",pch=0,xlab="Número de clusters",ylab="Silueta")
El método del elbow recomienda utilizar 3 clússteres para estos datos.
resultados <- rep(0, 7)
for (i in c(2,3,4,5,6,7))
{
fit <- kmeans(x, i)
resultados[i] <- fit$tot.withinss
}
plot(2:7,resultados[2:7],type="o",col="blue",pch=0,xlab="Número de clusters",ylab="tot.tot.withinss")
Si procedemos con la función kmeansruns bajos los criterios de la silueta media y Calinski-Harabasz obtenemos que para ambos casos el número óptimo de clústeres es de 3 y 2 respectivamente, aunque para Calinski-Harabasz tenemos que los valores que validan la opción de 3 clústeres es similar a la de 2.
fit_ch <- kmeansruns(x, krange = 1:10, criterion = "ch")
fit_asw <- kmeansruns(x, krange = 1:10, criterion = "asw")fit_ch$bestk## [1] 3fit_asw$bestk## [1] 2plot(1:10,fit_ch$crit,type="o",col="blue",pch=0,xlab="Número de clústers",ylab="Criterio Calinski-Harabasz")
plot(1:10,fit_asw$crit,type="o",col="blue",pch=0,xlab="Número de clústers",ylab="Criterio silueta media")
Procedemos pues a realizar el clusplot que claramente revela dos grupos y dos elipses, esta vez, en la dirección correcta.
data.kmeans.no_seasons.2_clusters <- kmeans(x, 2)
clusplot(x, data.kmeans.no_seasons.2_clusters$cluster, color=TRUE, shade=TRUE, labels=2, lines=0)
Al comparar la temperatura del aire media con otros factores, podemos ver una clara polarización de los datos cuando es contrastada con variables envueltas en las reglas de asociación, cosa que no ocurre con la precipitación.
par(mfrow=c(2,2))
plot(data.kmeans.no_seasons[c(1,3)], col=data.kmeans.no_seasons.2_clusters$cluster)
plot(data.kmeans.no_seasons[c(2,3)], col=data.kmeans.no_seasons.2_clusters$cluster)
plot(data.kmeans.no_seasons[c(4,3)], col=data.kmeans.no_seasons.2_clusters$cluster)
El model de tres clústeres se representa continuación.
data.kmeans.no_seasons.3_clusters <- kmeans(x, 3)
clusplot(x, data.kmeans.no_seasons.3_clusters$cluster, color=TRUE, shade=TRUE, labels=3, lines=0)
El contraste de variables se comporta similar de el dos clústeres.
par(mfrow=c(2,2))
plot(data.kmeans.no_seasons[c(1,3)], col=data.kmeans.no_seasons.3_clusters$cluster)
plot(data.kmeans.no_seasons[c(2,3)], col=data.kmeans.no_seasons.3_clusters$cluster)
plot(data.kmeans.no_seasons[c(4,3)], col=data.kmeans.no_seasons.3_clusters$cluster)
La precision es del 71.15%, pero en este caso para ambos modelos independienemete de los clústeres usados.
data.kmeans.no_seasons.2_clusters.table <- table(data.kmeans.no_seasons.2_clusters$cluster,data.kmeans.no_seasons$outcome)
data.kmeans.no_seasons.2_clusters.table##
## normal pandemic
## 1 477 37
## 2 559 383data.kmeans.no_seasons.2_clusters.table.precision = 100*(max(data.kmeans.no_seasons.2_clusters.table[1], data.kmeans.no_seasons.2_clusters.table[3])+ max(data.kmeans.no_seasons.2_clusters.table[2], data.kmeans.no_seasons.2_clusters.table[4]))/(sum(data.kmeans.no_seasons.2_clusters.table))
data.kmeans.no_seasons.2_clusters.table.precision## [1] 71.15385data.kmeans.no_seasons.3_clusters.table <- table(data.kmeans.no_seasons.3_clusters$cluster,data.kmeans.no_seasons$outcome)
data.kmeans.no_seasons.3_clusters.table##
## normal pandemic
## 1 468 37
## 2 280 270
## 3 288 113data.kmeans.no_seasons.3_clusters.table.precision = 100*(max(data.kmeans.no_seasons.3_clusters.table[1], data.kmeans.no_seasons.3_clusters.table[4])+ max(data.kmeans.no_seasons.3_clusters.table[2], data.kmeans.no_seasons.3_clusters.table[5]) + max(data.kmeans.no_seasons.3_clusters.table[3], data.kmeans.no_seasons.3_clusters.table[6]))/(sum(data.kmeans.no_seasons.3_clusters.table))
data.kmeans.no_seasons.3_clusters.table.precision## [1] 71.15385A pesar de que el modelo con los datos estacionales y cuatro clústeres denota un mejor rendimiento, consideramos que el de dos clústeres sin datos estacionales explica mejor el modelo. Esta conclusión la basamos en los siguientes hechos:
El modelo no estacional de dos clústeres se ajusta mejor a nuestro modelo de agrupación para los casos de dengue basado en la media.
El modelo no estacional hace uso de más variables para explicar las reglas mientras que el estacional delega la creación de reglas en la estación en la que se recogieron los datos.
Los datos estacionales son una constante independiente que sigue un patrón marcado y bien distribuido, cosa que no podemos afirmar del resto de las variables, que son medidas científicas de las condiciones climatológicas y morfológcas de las zonas estudiadas.
Los clusplot para los modelos estacionales invierte claramente la dirección de las zonas definidas y los clústeres de datos reales.
En resumen, creemos que la variable estacional explica el modelo desde una perspectiva poco científica y por ende será eliminada de los estudios que a continuación se producen.
Como método alternativo hemos seleccionado k-medoids basado en la distancia promedio de los puntos al punto real y central del cúmulo. Utilizaremos la función fviz_nbclust para realizar las estimaciones de los cúmulos. Calcularemos el número óptimos de clúster utilizando el método del codo y el de la distancia euclidiana. Como ya se ha mencionado anteriormente, omitiremos la variable estacional.
data.kmedoids <- data.kmeans.no_seasons
x <- data.kmedoids[,1:5]El número recomendado de clústeres es de 2
fviz_nbclust(x = x, FUNcluster = pam, method = "silhouette", k.max = 10, diss = dist(x, method = "euclidean"))
fviz_nbclust(x = x, FUNcluster = pam, method = "wss", k.max = 10, diss = dist(x, method = "euclidean"))
El gráfico de clusplot es similar al del ejercicio anterior con dos grupos de datos claramente señalados.
data.kmedoids.cluster <- pam(x = x, k = 2, metric = "euclidean")
clusplot(x, data.kmedoids.cluster$cluster, color=TRUE, shade=TRUE, labels=2, lines=0)
La precisión del model es del 71.15%, igual que para el método de kmeans.
data.kmedoids.cluster.table <- table(data.kmedoids.cluster$cluster,data.kmedoids$outcome)
data.kmedoids.cluster.table##
## normal pandemic
## 1 557 383
## 2 479 37data.kmedoids.cluster.table.precision = 100*(max(data.kmedoids.cluster.table[1], data.kmedoids.cluster.table[3])+ max(data.kmedoids.cluster.table[2], data.kmedoids.cluster.table[4]))/(sum(data.kmedoids.cluster.table))
data.kmedoids.cluster.table.precision## [1] 71.15385A continuación vamos a utilizar el metido Fuzzy clustering que hace uso del algoritmo c-means donde el centro del clúster es la media de la probabilidad de las observaciones de pertenecer a ese grupo.
data.kmedoids <- data.kmeans.no_seasons
x <- data.kmedoids[,1:5]Este método recomienda el uso de 3 clústeres.
data.fuz.output <- rep(0, 5)
for (i in c(2,3,4,5))
{
fit <- fanny(x = x, diss = FALSE, k = i, metric = "euclidean", stand = FALSE)
data.fuz.output[i] <- fit$coeff['normalized']
}
plot(2:5,data.fuz.output[2:5],type="o",col="blue",pch=0,xlab="Número de clusters",ylab="Coeficiente Dunn")
data.fuz.3_cluster <- fanny(x = x, diss = FALSE, k = 3, metric = "euclidean", stand = FALSE)
clusplot(x, data.fuz.3_cluster$clustering, color=TRUE, shade=TRUE, labels=2, lines=0)
Podemos observar que el rendimiento es de 71.49%.
data.fuz.3_cluster.table <- table(data.fuz.3_cluster$cluster,data.kmeans.no_seasons$outcome)
data.fuz.3_cluster.table##
## normal pandemic
## 1 302 122
## 2 256 261
## 3 478 37data.fuz.3_cluster.table.precision = 100*(max(data.fuz.3_cluster.table[1], data.fuz.3_cluster.table[4])+max(data.fuz.3_cluster.table[2], data.fuz.3_cluster.table[5])+max(data.fuz.3_cluster.table[3], data.fuz.3_cluster.table[6]))/(sum(data.fuz.3_cluster.table))
data.fuz.3_cluster.table.precision## [1] 71.49725En general, concluimos que los modelos de dos clústeres observados devengan una precisión del 71.15% en todos los métodos, mientras que el de 3 clústeres para el método de Fuzzy clustering nos aporta un 71.49% de precisión, algo superior que el obtenido para kmeans.
Vamos a crear ahora un modelo basado en árbol de decisiones, para ello utilizaremos el paquete CARET que nos permitirá optimizar los parámetros del modelo de forma recursiva.
data.tree <- data.kmeans.no_seasons
summary(data.tree)## ndvi precipitation air_temp_mean diur_temp_range
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.00000
## 1st Qu.:0.3239 1st Qu.:0.0255 1st Qu.:0.4005 1st Qu.:0.06621
## Median :0.4085 Median :0.0991 Median :0.5323 Median :0.10273
## Mean :0.4356 Mean :0.1172 Mean :0.5375 Mean :0.24174
## 3rd Qu.:0.5353 3rd Qu.:0.1793 3rd Qu.:0.6864 3rd Qu.:0.42235
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.00000
## humidity outcome
## Min. :0.0000 normal :1036
## 1st Qu.:0.4397 pandemic: 420
## Median :0.6116
## Mean :0.5752
## 3rd Qu.:0.7157
## Max. :1.0000Antes de proceder, creamos nuestro conjunto de datos de entrenamiento y pruebas.
set.seed(777)
train <- createDataPartition(data.tree$outcome, p = 0.66, list = FALSE)
data.tree.train <- data.tree[train,]
data.tree.test <- data.tree[-train,]
data.tree.train.x <- data.tree.train[,1:5]
data.tree.train.y <- data.tree.train[,6]
data.tree.test.x <- data.tree.test[,1:5]
data.tree.test.y <- data.tree.test[,6]data.tree.model.grid <- expand.grid(winnow = c(FALSE, TRUE), trials = seq(from = 1, to = 20, by = 2), model = c("rules"))
data.tree.model <- train(
outcome ~ .,
data = data.tree.train,
method = "C5.0",
tuneGrid = data.tree.model.grid,
metric = 'Accuracy'
)data.tree.model## C5.0
##
## 962 samples
## 5 predictor
## 2 classes: 'normal', 'pandemic'
##
## No pre-processing
## Resampling: Bootstrapped (25 reps)
## Summary of sample sizes: 962, 962, 962, 962, 962, 962, ...
## Resampling results across tuning parameters:
##
## winnow trials Accuracy Kappa
## FALSE 1 0.7196039 0.2612001
## FALSE 3 0.7128986 0.2730602
## FALSE 5 0.7154400 0.2747109
## FALSE 7 0.7152156 0.2705228
## FALSE 9 0.7156527 0.2825051
## FALSE 11 0.7170237 0.2830907
## FALSE 13 0.7142482 0.2884227
## FALSE 15 0.7151866 0.2892505
## FALSE 17 0.7150693 0.2895645
## FALSE 19 0.7148347 0.2886735
## TRUE 1 0.7271019 0.2854781
## TRUE 3 0.7185526 0.2882082
## TRUE 5 0.7237425 0.2930462
## TRUE 7 0.7252115 0.3020172
## TRUE 9 0.7252247 0.3054888
## TRUE 11 0.7246593 0.3034695
## TRUE 13 0.7245092 0.3094721
## TRUE 15 0.7251143 0.3091773
## TRUE 17 0.7257747 0.3116507
## TRUE 19 0.7254290 0.3092181
##
## Tuning parameter 'model' was held constant at a value of rules
## Accuracy was used to select the optimal model using the largest value.
## The final values used for the model were trials = 1, model = rules and winnow
## = TRUE.A continuación, extraemos el modelo más eficiente.
data.tree.model$bestTune## trials model winnow
## 11 1 rules TRUEAl siguiente gráfico nos revela la progresión de los modelos calculados según su precisión.
plot(data.tree.model)
Realizamos una predicción con los datos de pruebas, obsevando una precisión del 73.28%.
data.tree.predicted_model <- predict(data.tree.model, data.tree.test)
data.tree.predicted_model.confusionMatrix <- confusionMatrix(data.tree.predicted_model, data.tree.test$outcome)
data.tree.predicted_model.confusionMatrix## Confusion Matrix and Statistics
##
## Reference
## Prediction normal pandemic
## normal 309 89
## pandemic 43 53
##
## Accuracy : 0.7328
## 95% CI : (0.6914, 0.7713)
## No Information Rate : 0.7126
## P-Value [Acc > NIR] : 0.1727
##
## Kappa : 0.2779
##
## Mcnemar's Test P-Value : 8.975e-05
##
## Sensitivity : 0.8778
## Specificity : 0.3732
## Pos Pred Value : 0.7764
## Neg Pred Value : 0.5521
## Prevalence : 0.7126
## Detection Rate : 0.6255
## Detection Prevalence : 0.8057
## Balanced Accuracy : 0.6255
##
## 'Positive' Class : normal
## Podemos decir que este modelo realiza un mejor trabajo realizando predicciones del tipo normal en contraste con las pandemicas. Esto es un reflejo de la composición de los datos y su clasificación, donde existen más valores que dan como resultado un diagnóstico normal. Esto puedo verse reflejado en el alto nivel de sensibilidad ( true positive ) del test en comparación con la especificidad ( true negative ).
A continuación, vamos a efectuar un análisis PCA para descartar variables que puedan estar aportando poco al modelo final. Este análisis nos permitirá simplificar el modelo con un bajo coste en la precisión final del modelo aunque con una ganancia en rendimiento a la hora de calcular el modelo a utilizar.
data.tree.pca.analysis <- prcomp(data.tree[,1:5], center = TRUE, scale. = TRUE)
summary(data.tree.pca.analysis)## Importance of components:
## PC1 PC2 PC3 PC4 PC5
## Standard deviation 1.3922 1.2904 0.8815 0.60091 0.50865
## Proportion of Variance 0.3876 0.3330 0.1554 0.07222 0.05174
## Cumulative Proportion 0.3876 0.7206 0.8760 0.94826 1.00000Si observamos al resultado del modelo, podemos deducir que son los dos primeros componentes quienes explican mejor el modelo. De la rotación podemos observar que ndvi, air_temp_mean y diur_temp_range tienen más peso a la hora de calcular el modelo.
data.tree.pca.analysis$rotation## PC1 PC2 PC3 PC4 PC5
## ndvi 0.61949453 -0.09773145 0.3302020 -0.5119503 0.4853335
## precipitation 0.23064993 -0.52108003 -0.7144128 0.2330397 0.3325403
## air_temp_mean -0.41994894 -0.45286855 0.5264852 0.3566120 0.4628112
## diur_temp_range 0.61929755 -0.07972356 0.3026169 0.6311647 -0.3466537
## humidity -0.05604084 -0.71237518 0.1087080 -0.3975655 -0.5652480Si realizamos el plot de los dos primeros componentes, podemos observar que es PC1 quien explica mejor el modelo y que existe una mayor descentralización de los dos grupos. Además, podemos observar como ndvi_, air_temp_mean y diur_temp_range son las encargadas de polarizar el diagrama desde las esquinas éste.
ggbiplot(data.tree.pca.analysis,ellipse=TRUE, obs.scale = 1, var.scale = 1, groups=data.tree$outcome, alpha = 0)
ggbiplot(data.tree.pca.analysis,ellipse=TRUE, choices = c(2,3), obs.scale = 1, var.scale = 1, groups=data.tree$outcome, alpha = 0)
Vamos ahora a repetir el modelo de predicción basado en el árbol de decisiones del paquete CARET.
data.tree.pca <- data.tree
summary(data.tree.pca)## ndvi precipitation air_temp_mean diur_temp_range
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.00000
## 1st Qu.:0.3239 1st Qu.:0.0255 1st Qu.:0.4005 1st Qu.:0.06621
## Median :0.4085 Median :0.0991 Median :0.5323 Median :0.10273
## Mean :0.4356 Mean :0.1172 Mean :0.5375 Mean :0.24174
## 3rd Qu.:0.5353 3rd Qu.:0.1793 3rd Qu.:0.6864 3rd Qu.:0.42235
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.00000
## humidity outcome
## Min. :0.0000 normal :1036
## 1st Qu.:0.4397 pandemic: 420
## Median :0.6116
## Mean :0.5752
## 3rd Qu.:0.7157
## Max. :1.0000set.seed(777)
train <- createDataPartition(data.tree.pca$outcome, p = 0.66, list = FALSE)
data.tree.pca.train <- data.tree.pca[train,]
data.tree.pca.test <- data.tree.pca[-train,]
data.tree.pca.train.x <- data.tree.pca.train[,1:5]
data.tree.pca.train.y <- data.tree.pca.train[,6]
data.tree.pca.test.x <- data.tree.pca.test[,1:5]
data.tree.pca.test.y <- data.tree.pca.test[,6]Entrenamos un modelo con CARET indicando que queremos un preprocesamiento de los datos con PCA. CARET se encargará de seleccionar las variables óptimas para el modelo con variables reducidas.
data.tree.pca.model.grid <- expand.grid(winnow = c(FALSE, TRUE), trials = seq(from = 1, to = 15, by = 2), model = c("rules"))
data.tree.pca.model <- train(
outcome ~ .,
data = data.tree.pca.train,
method = "C5.0",
tuneGrid = data.tree.pca.model.grid,
metric = 'Accuracy',
preProcess = c('pca')
)data.tree.pca.model## C5.0
##
## 962 samples
## 5 predictor
## 2 classes: 'normal', 'pandemic'
##
## Pre-processing: principal component signal extraction (5), centered (5),
## scaled (5)
## Resampling: Bootstrapped (25 reps)
## Summary of sample sizes: 962, 962, 962, 962, 962, 962, ...
## Resampling results across tuning parameters:
##
## winnow trials Accuracy Kappa
## FALSE 1 0.7067199 0.2181694
## FALSE 3 0.7066854 0.1989197
## FALSE 5 0.7092433 0.2144911
## FALSE 7 0.7111387 0.2258201
## FALSE 9 0.7105611 0.2289777
## FALSE 11 0.7092170 0.2272623
## FALSE 13 0.7094411 0.2278809
## FALSE 15 0.7092170 0.2272623
## TRUE 1 0.7098970 0.2082943
## TRUE 3 0.7081524 0.2150357
## TRUE 5 0.7102417 0.2178224
## TRUE 7 0.7099807 0.2233506
## TRUE 9 0.7115944 0.2348943
## TRUE 11 0.7114823 0.2348682
## TRUE 13 0.7117064 0.2354868
## TRUE 15 0.7114823 0.2348682
##
## Tuning parameter 'model' was held constant at a value of rules
## Accuracy was used to select the optimal model using the largest value.
## The final values used for the model were trials = 13, model = rules and
## winnow = TRUE.data.tree.pca.model$bestTune## trials model winnow
## 15 13 rules TRUEplot(data.tree.pca.model)
Finalmente vemos que la precisión del modelo es de 71.46%, un 1.82% inferior que el anterior sin PCA.
data.tree.pca.predicted_model <- predict(data.tree.pca.model, data.tree.pca.test)
data.tree.pca.predicted_model.confusionMatrix <- confusionMatrix(data.tree.pca.predicted_model, data.tree.pca.test$outcome)
data.tree.pca.predicted_model.confusionMatrix## Confusion Matrix and Statistics
##
## Reference
## Prediction normal pandemic
## normal 316 105
## pandemic 36 37
##
## Accuracy : 0.7146
## 95% CI : (0.6725, 0.754)
## No Information Rate : 0.7126
## P-Value [Acc > NIR] : 0.483
##
## Kappa : 0.1851
##
## Mcnemar's Test P-Value : 1.024e-08
##
## Sensitivity : 0.8977
## Specificity : 0.2606
## Pos Pred Value : 0.7506
## Neg Pred Value : 0.5068
## Prevalence : 0.7126
## Detection Rate : 0.6397
## Detection Prevalence : 0.8522
## Balanced Accuracy : 0.5791
##
## 'Positive' Class : normal
## Este modelo ha aumentado su capacidad para realizar predicciones del tipo normal, pero ha perdido capacidad en lo referente a predicciones del tipo pandemic. Esto ha ocurrido como consecuencia de la eliminación de variables que, aunque aportaban poco al modelo, favorecían a la predicción en favor de la clase negativa o pandemic.
Como bien se ha mencionado, existe un decremento en la precisión del modelo con PCA del 1.82%. Este es debido a que estamos trabajando con menos variables y por consecuente el modelo será menos preciso. Sin embargo, hay que tener claro que el propósito de PCA no es el de mejorar el rendimiento del modelo sino el de simplificarlo con el menor coste posible, que es lo que se ha logrado con este procedimiento.